REPORT  79-27 

Some  Bessel  function  identities  arising 
in  ice  mechanics  problems 


o 

00 

o 

Q 


51  ho(/» 


m  _ 

eT  \fx2  +  y 


0 


JSx  +  iny  Hw  —  - — - 

e  dxdy  -  +  +  pi 


4/ 


1__  e~i*i-irn  d|dr,  =  rtW<"(/?e^  x/F+T5) 


42  +  rj2  +  /J 


This  doc  ument  ia<e»~ 
lor  public  r«*lco»o  cn  -  . 
■  distribution  ie  tmliauitcc; 


iporr  v<vr 


79  12  31  034 


CRREL  Rejsoft-79-27 


$ome  Bessel  function  identities  arising 
in  ice  mechanics  problems  * 


)  ^hunsukeyTakagi 


Prepared  for 

DIRECTORATE  OF  MILITARY  PROGRAMS 
OFFICE,  CHIEF  OF  ENGINEERS 


UNITED  STATES  ARMY 
CORPS  OF  ENGINEERS 

COLD  RECIONS  RESEARCH  AND  ENGINEERING  LABORATORY 
HANOVER,  NEW  HAMPSHIRE,  U  S  A 


Approved  for  public  release,  distribution  unlimited 


_ Unclassified _ 

SECURITY  CLASSIFICATION  OF  THIS  PACE  (Whan  Data  E ntarad) 


REPORT  DOCUMENTATION  PAGE 

READ  INSTRUCTIONS 

BEFORE  COMPLETING  FORM 

t.  REPORT  NUMBER  2.  GOVT  ACCESSION  NO. 

CRREL  Report  79-27 

3.  RECIPIENT’S  CATALOG  NUMBER 

4.  TITLE  rand  Submit) 

SOME  BESSEL  FUNCTION  IDENTITIES  ARISING  IN  ICE 
MECHANICS  PROBLEMS 

5.  type  of  report  a  period  COVERED 

6.  PERFORMING  ORG.  REPORT  NUMBER 

7.  author^; 

Shunsuke  Takagi 

8.  CONTRACT  OR  GRANT  NUMBER^ 

9.  performing  organization  name  and  aooress 

U.S.  Army  Cold  Regions  Research  and  Engineering  Laboratory  ^ 
Hanover,  New  Hampshire  03755 

10.  PROGRAM  ELEMENT.  PROJECT.  TASK 
AREA  A  WORK  UNIT  NUMBERS 

DA  Project  4A1 61 102AT24 

Scientific  Area  A,  Work  Unit  005 

11.  CONTROLLING  OFFICE  NAME  AND  ADDRESS 

Directorate  of  Military  Programs 

Office,  Chief  of  Engineers 

Washington,  D.C.  20314 

12.  REPORT  DATE 

November  1979 

13.  NUMBER  OF  PAGES 

14 

1*.  MONITORING  AGENCY  NAME  4  ADDRESS (It  dlllarant  front  Controlling  Otllco) 

15.  SECURITY  CLASS,  (of  thla  report) 

Unclassified 

15a.  DECLASSIFICATION/ DOWNGRADING 
SCHEDULE 

16.  DISTRIBUTION  STATEMENT  (of  thla  Report) 

Approved  for  public  release;  distribution  unlimited. 

D  D  C 

npaEIME. 

JAN  2 

n 

17.  DISTRIBUTION  STATEMENT  (ol  the  obatract  entered  In  Block  30,  II  dlllarant  from  Report) 

UUEESIETOTb 

E 

I 

• 

16.  supplementary  notes 

19.  KEY  WORDS  (Continue  on  r«v«r*«  aid*  If  n*c*aa*ry  and  Identify  by  block  number) 

Bessel  functions 

Fourier  integrals 
|  Ice  mechanics 

Winkler  foundations 

,  \ 

VU  ABSTRACT  fCantBaue  an  rev  area  aid*  ft  ntcttMqr  and  Identity  by  block  number) 

Some  Bessel  function  identities  found  by  solving  problems  of  the  deflection  of  a  floating  ice  plate  by  two  different 
methods  are  rigorously  proved.  The  master  formulas  from  which  all  the  identities  are  derived  are  in  a  Fourier  recip¬ 
rocal  relationship,  connecting  a  Hankel  function  to  an  exponential  function.  Many  new  formulas  can  be  derived  from 
the  master  formulas.  The  analytical  method  presented  here  now  opens  the  way  to  study  a  hitherto  impossible  type 
of  problem -the  deflection  of  floating  elastic  plates  of  various  shapes  and  boundary  conditions. 

DO  %'%?n  1473  EDITION  OF  I  NOV  CS  IS  OBSOLETE 


Unclassified 


SECURITY  CLASSIFICATION  OF  THIS  PAGE  Dal*  E ntarad) 


i  a---* 


PREFACE 

This  report  was  prepared  by  Dr.  Shunsuke  Takagi,  Research  Physical  Scientist,  of  the  Physical 
Sciences  Branch,  Research  Division,  U.S.  Army  Cold  Regions  Research  and  Engineering  Laboratory 
Funding  for  this  research  was  provided  by  DA  Project  4A161 1 02AT24,  Research  in  Snow,  Ice  and 
Frozen  Ground,  Scientific  Area  A,  Properties  of  Cold  Regions  Materials,  Work  Unit  005,  Thermo- 
physics  and  Mechanics  of  Cold  Regions  Materials. 


SOME  BESSEL  FUNCTION  IDENTITIES 
ARISING  IN  ICE  MECHANICS 
PROBLEMS 

Shunsuke  Takagi 


INTRODUCTION 

By  solving  several  problems  of  the  deflection  of  a  floating  ice  plate  (i.e.  the  deflection  of  a  plate 
on  a  continuous  elastic  foundation  as  formulated  by  Winkler  1 867)  by  two  different  methods,  Kerr 
(1977a,  1978)  presented  a  number  of  equality  relationships  among  Bessel  functions.  In  this  report, 
the  analytical  derivation  of  his  formulas  is  presented.  His  formulas  (six  in  all)  reduce  to  the  following 
two  master  formulas  expressed  in  the  Fourier  reciprocal  relationship: 


H^^e  2\^+ 


\/a 2  +  £2  jd£  = 


< 2  +  02  -  i(x  dx 
s/x2  +B2 


Hn>U2#TF, 


where  hJj'  ^  ( )  is  the  Hankel  function  of  zeroth  order;  x  and  £  are  real  and  a  nonnegative  such  that 
a2+£2  0;  and  P  is  complex  such  that  (3  =£  0,  larg/?l  <  ir/2,x2+p2  ^  0.  By \Jz  we  mean  such 

a  branch  as  R e\/z  2  0.  The  master  formulas  may  be  transformed  to  double  integrals  that  are  sym¬ 
metric  with  regard  to  arguments: 
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In  the  following,  we  prove  the  above  formulas,  transform  the  two  master  formulas  to  derive  all 
the  formulas  introduced  by  Kerr  (1977a,  1978)  as  well  as  some  new  ones,  and  finally  show  that  the 
analytical  method  presented  here  lays  the  foundation  for  building  mathematical  machinery  for 
solving  problems  for  various  shapes  of  floating  ice  plates  under  various  boundary  conditions. 
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PROOF  OF  FORMULA  (1) 


Use  of  Barnes’  integral  representation  of  H^1  *  (z), 

■M-ff*'  (5) 

where  e  is  any  positive  number  and  larg  (-/Z)l  <  tt/2  (Watson  1962,  p.  192),  transforms  the  single 
integral  on  the  left-hand  side  of  eq  1, 

I,  =  ni  jT  eixt  H*,' 
to  the  repeated  integrals 

•’  -  s-L  r^(  i  'P7*T  <7) 

The  absolute  convergence  of  the  integral  (5)  carries  over  to  (7),  because  the  condition  larg(-/z)l 
<  tj/2  in  (5)  transforms  to  larg/3l  <  ir/2  in  (7),  which  is  one  of  the  prerequisites  in  the  master  for¬ 
mulas.  Therefore,  the  order  of  integration  in  (7)  may  be  exchanged.  Moreover,  restricting  the  orig¬ 
inal  range  of  a,  which  is  a  >  0,  too  >  0,  we  let  |  =077  in  (7)  so  that  it  becomes 


To  evaluate  the  internal  single  integral  in  (8) 
1,  =  J  e'ax/  (1  +z2)s  dz 

»  —00 


.2  vs 
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M, 


by  the  contour  integral  method,  we  first  note  that  the  original  range  of  x  (i.e.  -  00  <  x  <  °°)  may 
be  restricted  to  0  <  x  (where  the  case  x  =  0  is  excluded),  because  I]  is  an  even  function  of  x,  as 
the  right-hand  side  of  (6)  shows.  Consider  the  contour  in  Figure  1  that  starts  at  origin  0,  goes  along 
the  positive  real  axis  to  A  (i.e.  z  =  °°),  takes  a  90°  turn  along  the  infinitely  large  circle  to  reach  B 
( i.e.  z  =  /°°),  comes  down  along  the  imaginary  axis  to  C  (i.e.  z  =  /),  makes  a  360°turn  along  an 
infinitely  small  circle  clockwise  around  C,  goes  upward  along  the  imaginary  axis  to  reach  D  (i.e. 
z  =  /o°),  takes  a  90°  turn  along  the  infinitely  large  circle  to  reach  E  (i.e.  z  =  -00),  and  finally  reaches 
origin  0,  thus  completing  a  circuit.  No  singularity  of  the  integrand  exp (iaxz)  (1  +z2)s  exists  inside 
this  closed  contour.  Among  the  integrals  along  the  paths  mentioned  above,  the  integrals  along  AB 
and  DE  are  equal  to  zero,  provided  that  a  >  0.  The  integral  around  C  is  also  equal  to  zero.  There¬ 
fore,  on  the  condition  that  a  >  Oandx  >  0,  we  have 


M,  = 


e/o*z  (1  +Z2)J  d/f 


where  Cj  and  C2  are  the  initial  and  terminal  points  of  the  infinitely  small  circle  around  C.  Letting 
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Figure  1.  Transformation  of  integral  M1 
on  the  i-piane. 


z  -  it,  where  t  is  real,  M,  reduces  to 


M,  =  /(I  -e-2w/s)  J  e~axt  (1  -f2)5  dr. 

We  now  let  5  be 
s  =  -V$  +  /'p; 

i.e.  letc  in  (5)  be  A,  where  p  is  a  real  number,  in  order  to  integrate  by  use  of  the  formula  in 
Watson  (1962,  p.  172), 


K,W  =  f  d t, 

r(t» + 'A)  J  i 


(9) 


which  is  valid  when  Re(K  +  '/2)  >  0  and  largzl  <  n/2.  These  two  conditions  are  satisfied  when  we 
\etv-'A  =  sand  z  =  ax  to  integrate  M,.  Thus,  letting 

(1  -r2)5  =  e"/s  (f 2  - 1)1, 

M,  integrates  to 

M,  =  -2  sin  ns — CIlLLLKj  +  «  (ax). 

In  this  way,  (7)  transforms  to  a  single  integral: 

-}A  +  oof 
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•  /*  r2(-s)/^.\2  sin  ns — P(5  +  1).  K  »(ax)  ds.  V  vo 


Changing  the  Gamma  function  of  the  negative  argument  to  the  positive  argument  by  the  reflection 
formula, 


T(-s)  = 


T(1  +s)  sin  ns  ' 
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Ij  becomes 


p-Yi  +  °°  / 

1  =-1^h|  * 

1  w'v  x  J -Yt  -«/  sin  us  r(s  +  1) 
Taking  the  residues,  Ij  integrates  to 

Replacing  K„  +  yl(ax)  with 


i.e.,  a  formula  found  in  (8,  p.  1 83] ,  I,  becomes 


The  order  of  the  summation  and  integration  may  be  exchanged,  and  we  get 


I 


I 


dr 

4t  v 


Letting  t  =  £'2  ,  this  becomes 


00 

I,  =  2 afifo  exp  (-r2  -  g.V+.g-2Lg2)  dl  (1 0) 

To  integrate  (10),  we  introduce  a  lemma: 

f  e^2«2-*‘2«  =  e'2*1^, 

Jo  2>i 

provided  that 

m r-  —  <  argp  <  r»ff+- 
4  4 

where  n  is  an  integer. 

When  n  is  in  the  above  range,  the  above  integral  is  convergent.  To  prove  the  lemma,  we  first  note 
that  the  integral 


transforms  to  L,  =  e_2M  N,, 


where  N,  =  e~^u*  *  d£. 


Letting  £  =  1  /(ait/) 

and  changing  the  resulting  contour  0  ~  /a-1 00  to  0  ~  °°,  we  get 


N ,  =  f  e _*L. 
Jo  ,,„2 


Addition  of  (1 2)  and  (13)  yields 


Letting 


^i-r1  =  t, 


u  00 

*-*jL 


e“r  d t. 


Changing  the  range  of  integration  to  the  one  from -°°  to +°°,  Nj  integrates  to  Nj  =  n/(2u).  Sub¬ 

stituting  this  value  into  (1 1 ),  the  lemma  is  proved. 

Letting  /a  be 

H  = 

in  the  lemma,  (10)  is  integrated,  because  u  above  is  obviously  in  the  range  prescribed  before.  Thus, 
under  the  conditions  x  *=  Oando  0,  Formula  (1)  is  proved.  Applying  the  analytical  continua¬ 
tion,  the  condition  a  *  0  is  extended  to  the  condition  <72  +  02  ¥=  0.  Because  the  integral  is  con¬ 
vergent  at  x  =  0,  the  condition  x  =£  0  may  be  removed.  The  proof  is  thus  completed. 


PROOF  OF  FORMULA  (2) 

Although  Formula  (2)  is  the  Fourier  inverse  of  Formula  (1),  we  show  an  independent  proof  in 
view  of  the  importance  of  the  formula. 

On  the  assumption  that  <72  +  £2  +  Oand0  =£  0,  letting 

a  =  rcosa 

{  =  rsina  (1 

x  =  /JSinhz 
-  —  <  a  <  — , 


c 


the  integral 


where  the  contour  L  is  a  curve  on  the  complex  z  -  u  +iv  plane  (Fig.  2)  defined  by 
/  =  Arcsinh  (*//?) 


(15) 


(16) 


with  parameter  x  in  the  range  of  <  x  < 

Letting x  =  0,  we  have/  =  0.  Therefore,  contour  L  passes  through  the  origin.  When  x  -»  +  » 
or  -°°,z  asymptotically  approaches  log  (2x)  -  log  0,  or  -  log  (-2x)  +  Iog0,  respectively;  in  other  words, 
the  imaginary  part  v  of  complex  variable  z  satisfies  the  conditions 
lim  v  +  arg/3  =  0 
*-*•+«> 

(18) 

lim  y-arg/J  =  0 


The  curve  L+  defined  for  the  case  0  <  arg0  <  jt/2  is  shown  in  Figure  2.  The  curve  L_  defined  for 
the  case  -tt/2  <  arg/3  <  0  is  symmetrical  with  L+  with  regard  to  the  real  axis. 

Letting  u  be  the  real  part  of  z,  we  find  that,  as  x  -►  °°  or  -  »,  the  real  part  of  -r  /3Cosh  (z  +  /<*) 
approaches  asymptotically  either 


1_ 

2 


r l/Jle“cos  (arg^  +  r'  +  a) 


or 


-1  rl^le""cos  (arg/3-v-a), 

respectively.  Because  the  power  of  the  exponent  in  (16)  must  remain  negative  as  Ixl  -*■  °°,  the 
conditions 


7T 

2 


< 


arg/3  +  lim  v  +  a  <  — 

x-»+~ 


-  -  <  arg)3-lim  v-a  <  — 

2  2 
X~*  -<*> 

must  be  satisfied.  The  conditions  are  satisfied  by  (14)4  and  (1 8). 


(19) 
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Figure  2.  Curve  L+  defined  for  the  case  0  <  argl 3  <  v/2  in  the  com¬ 
plex  i  =  \}  +  \v  plane  is  shown  for  parameter  x  ranging  over-  °°  <  x 

< 


Introducing  t  defined  by 

t  =  z+ia  (20) 

l2  becomes 


■j f 

Jl  +  la 


exp  (-rp  Cosh  t)  d t. 


The  contour  L  +/'o  may  be  moved  to  the  real  axis,  because  doing  this  is  tantamount  to  letting 

which  is  evidently  permissible  in  the  prescribed  range  of  0  shown  below  (2).  Thus 

we  have 


TRANSFORMATIONS  OF  FORMULAS  (1)  AND  (2) 

We  list  the  formulas  found  in  Erdelyi  et  al.  (1954)  and  Kerr  (1978, 1977b)  that  can  be  derived 
from  Formulas  (1 )  and  (2).  The  formulas  we  use  in  the  following  transformations  are 


7 


which  is  found  by  exchanging  x  and  a  in  Formula  (1 ),  and 


£  e-ovT2  +  02  -ibx  —4* - =  rr/H*,1  *  (/3  e  2  Ja2  +  62), 

which  is  found  by  changing  £  to  b  in  Formula  (2). 

Assuming  0  to  be  a  positive  number,  (23)  transforms  to  a  real  integral : 


(24) 


fo  K0  (M7^ 


f)  cos  (o|)  d£  = 


* _ P-V^7+  02* 


2^W2 


(25) 


Letting  a  -  0,  (25)  may  become 

I  K0  (0x^77)  di  =  ie-»\ 

J-oo  ft 


(26) 


Letting  £  -  r)-y,  and  expressly  specifying  that  a  may  be  positive  or  negative,  (23)  may  be  re¬ 
written  as 


£  et/aV  H("(0e2V^T ^j7) 


drj  =  2  t±/ay-\/a2  +  e2x 

i'Ja2  +02 

We  transform  this  to  several  forms.  When  0  is  real,  (27)  becomes 


(27) 


£  elto’>  K0  ^0v42 +(y-rj)2j 


dr? 


_£ _ -t/ay  -Va2  +  p2x 


v£2+02 


(28) 


Letting  0  =  b  exp  (ff//4)  with  the  restriction  b  >  0,  (27)  becomes 


£  etlar>  Jker  (bVx2  +  (y-tj)2)  +/kei  (by/x2  +  (y  -r?)2j  Jdr?  =  n  '62*  *  '< 

V^2  +/fc2 


(29) 


Letting 

w2+/62  =  p+/<7  (30) 

and  adding  and  subtracting  the  plus  expression  and  the  minus  expression  in  (29),  we  find  two  integrals: 


I3  = - 1 _ f  elvy  -  'y'y/i2+t2  dy 


Use  of  the  integral 


Jvy-  'ytJr  +  e2 , 


J2+r?2+/32 


reduces  l3  to  the  right-hand  side  of  Formula  (3).  To  prove  (35),  note  that 


j  ehy  -  l y  'y/t2  *H2  dy  =  2  J £  e->\45  +  g5  , 


cos  r)y  6y. 


Although  Formula  (4)  is  the  Fourier  inverse  of  Formula  (3),  we  show  an  independent  proof  i 
view  of  the  importance  of  the  formula.  The  rewritten  Formula  (2), 


yjr  *(r 


reduces  to  Formula  (4)  by  the  substitution  of  the  Fourier  inverse  of  (35): 


■lyiyA2*02  =  \k2+&2 


*  J~t2W+fi2 


2 - e'lyn  dr?. 


Formula  (36)  may  otherwise  be  proved  by  showing  that 


e"y/  dz  =  le-*1'1 


2  2 
i  +a* 


by  use  of  the  contour  integral  method,  where  Re  a  =£  0. 


ADDITIONAL  DERIVATION  FROM  THE  MASTER  FORMULAS 

We  prove  the  formula 


X.—-LrrHnJ)i(  fr2  J,1  +  |2  ^(fle2  v42  +  (y-j;)2  )df  =  ~  H^/Wt"**)2  +  v2), 

Vo2+i2  ^ 

(38) 
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Differentiating  (24)  with  regard  too,  we  find 


-/**  dx  =  n(Sa  Ho 


’  >■{  (Je2v^2). 


Substitution  of  (3ew/2  ,/a2  +  ^  from  (40)  transforms  the  single  integral, 

'• '  L  7^= «[,;),( )<■ 
va  +r 

to  the  reperted  integrals, 

>4  ^2+(y-l?  )d|  jT  t*^*?-*'  dt, 

which,  on  changing  the  order  of  the  integration,  becomes 


■  4  =  I  e 


uw2  +  02 


d' ^  e/tf  H*,1 ) (/?e  2 \/x2  +  (y-$2) )d{. 


Letting  {  -  rj-y,  and  changing  a  to  -f,  (23)  becomes 
f  e~'tT>  H|,1)(/3e2\42  +  (y-T?)2  )dt?  = - ? — ( 


-/f>*  -V^2  ♦  fl2  X 


Using  the  last  integral  to  carry  out  the  internal  integration,  l4  becomes 


J _ e-(»  *  x)  Vt 2  +  (T  -ity  , 


which  integrates  to 


U  =  27T  Ht,1  >( 


^e2V(o+Ar)2+y2' 


by  use  of  Formula  (2).  Combining  (41 )  and  (42),  (38)  is  proved. 
Letting  0  =  1  in  (38)  we  find 


DEFLECTION  OF  THE  FLOATING  ICE  PLATE 


Expressed  in  nondimensional  form,  the  differential  equation  governing  the  deflection  w  of  a 
floating  ice  plate  sustaining  a  concentrated  load  P  at  the  origin  (x  =  0 ,y  -  0)  is 

V4w  +  R'  =  P6(x)6(y),  (44) 

whereV2  is  the  Laplacian  operator  and  5( )  the  delta  function.  The  solution  of  (44)  for  an  infinite 
plate  is 


w{*,y)  =  -  ^-kei\42+>'2,  (45) 

as  was  shown  by  Wymann  (1950)  by  examining  the  nature  of  the  solution  of  the  homogeneous  form 
of  (44).  His  proof,  however,  does  not  exactly  show  that  (45)  is  the  solution  of  the  inhomogeneous 
equation  (44).  We  now  can  prove  this  by  direct  substitution  of  (45)  into  (44). 

Use  of  Formula  (4)  enables  us  to  derive 


kei\X2  +y2  = 


4-ff  - ! - e+lx*  *  dgdti. 

2nJ~  J—  (|2+t?2)2  +  1 


(46) 


With  the  use  of  (46),  substitution  of  (45)  reduces  the  left-hand  side  of  (44)  to 


e+ixt  v  iy-q  d{dj?i 


which  by  the  property  of  the  delta  function  (Gel'fand  and  Shilov  1964)  is  the  right-hand  side  of  (44). 

Use  of  the  reciprocal  Fourier  relationships  developed  in  this  report  opens  the  possibility  of  building 
analytical  machinery  for  solving  the  generic  equation 


V4w  +  w  =  P 6(x-x0)  5(y-y0) 


(47) 


for  the  deflection  of  floating  plates  of  various  shapes  and  of  various  boundary  conditions,  which  has 
hitherto  been  impossible.  (Currently  only  the  image  method  is  used;  see  Kerr  (1963  or  1977b).  It 
is  especially  encouraging  to  note  that  the  solution  w(x,y)  of  (47)  is  a  generalized  function  (Gel'fand 
and  Shilov  1964). 
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